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by 


J. R. Herman and S, Chandra 
Laboratory for Space Sciences 
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Greenbelt, Maryland 

ABSTRACT 

The electron continuity equation and the heat conduction 
equations for electrons, ions, and neutral species are solved 
simultaneously for specified boundary conditions and values 
of the solar flux. These solutions have yielded a set of self- 
consistent steady state profiles for a variety of midday solar 
conditions that are in agreement with observational data. To 
produce these profiles only the boundary conditions and solar 
flux are used as variable parameters of the problem. It is 
shown that the neutral temperature and the resulting neutral 
gas composition pla,y a dominant role in determining the charged 
particle density and temperature profiles. This leads to a 
picture of the solar cycle variation where composition changes 
in the neutral gas (introduced at the lower boundary) must toe 
combined with the solar flux variation to produce physically 
reasonable results. Of the five neutral gases used, (N„ 0 o , 

0, He, H) the atomic oxygen density variations are shown to be 
most effective in producing the observed trends for all the 
temperatures and densities in the E and F regions, 

*To be published in PLANETARY AND SPACE SCIENCE 
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INTRODUCTION 

The importance of obtaining a simultaneous solution of the 
equations of continuity, momentum, and heat transport for the 
electron, ion, and neutral gases has been recognized for some 
time. The processes influencing the general behavior of the 
E and F regions arising from changes in solar activity have 
a strong interdependence. A change in solar EUV radiation 
directly alters the photoionization rate and photoelectron 
energy spectrum. These in turn produce a sequence of secondary 
or indirect effects whose influence on the ionosphere is quite 
complex. For example, a fraction of the solar EUV radiation 
is absorbed by the neutral atmosphere and transformed into 
internal energy. The release of this energy by super-elastic 
collisional heating results in composition and temperature 
changes of the ionizable constituents; thus indirectly effecting 
the photoionization rate. The charged particle loss and 
transport processes are similarly effected by a series of cyclic 
chain reactions involving density, composition, and temperature 
of all the ionospheric constituents. Figure I shows a rough 
schematic illustrating the basic coupling between the various 
components making up the ionospheric system. In solving the 
equation of continuity or heat transport it has been usual to 
assume a suitable model for the unknown quantities appearing 
in the equations. This approach, while quite useful in under- 
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standing the general physical processes governing the E and F 
regions, masks the essential processes needed for studying 
the behavior due to changes in solar activity . 

The recent largely successful effort to make simultaneous 
measurements of the temperatures and densities of the ionizable 
and neutral constituents in the E and F regions (Geoprobe) , 
has increased the importance of obtaining a corresponding 
simultaneous theoretical description for the upper atmosphere 0 
It is hoped that a clearer understanding of the cause and effect 
relationships can be obtained through such an effort. 

In this work a simultaneous solution is described of the 
continuity and heat conduction equations for an ionospheric 
plasma consisting of electrons, atomic oxygen ions, and the 
neutral gas species N 2 9 0 2 , 0, He, and H , The solutions are 
presented for steady state conditions obtained by integrating 
the time dependent differential equations with time independent 
driving functions (stationary sun) until the transients introduced 
by the arbitrary choice of initial conditions vanish. 

Various studies have been conducted based on these equations 
[see Rishbeth (1968) , Yonezawa (1966) , and Evans (1966) for general 
references]? using slightly modified forms of the transport 
coefficients described in the earlier works of Chapman (1952) , 
Nicolet (1961) ,, and Spitzer (1965), More recently Dalgarno 
et . al. (1963) and Banks (1966) have listed many of the detailed 
expressions for the transport coefficients and driving functions 
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used in the ionosphere. Even though each of these works essentially 
analysed a single equation, the results demonstrate that the 
equations are sufficient to account for many features of the 
ionospheric plasma. 

In this paper we have followed the forms of the equations 
as given by Harris and Priester (1962) for the neutral temperature, 
Shimizaki (1965) for the electron continuity equation, and 
Banks (1966) for the charged particle temperatures. Because of 
these earlier works only an outline of the equations is presented 
instead of the detailed derivations. 

ELECTRON TEMPERATURE 

In order to describe the electron temperature profiles 
for various solar conditions it is necessary to consider the 
dominant effects of heat conduction and net local heat production. 
These effects correctly give the magnitude and shape of the 
temperature profile throughout the region of interest for 
midday conditions. Other terms can be included in the heat 
conduction equation to further refine the model, such as the 
effect of convection and non-local heating. However, such 
refinements have little effect on the steady state solutions 
and their interpretations as presented here. On the other hand, 
many of the interesting features of the diurnal behavior of the 
ionosphere arise from the inclusion of these extra terms. 

Accordingly, the following form of the heat conduction 
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equation is used 


for the electron temperature, 
bT 

n e kC v JT - [K e y = V L 


e 


( 1 ) 


If the main heat flow is along the earth’s magnetic field lines, 
and if the flat earth approximation is used, then in terms of 
the altitude z, the time t, and the magnetic dip angle I, the 

electron temperature profile T (z) is given by [Geisler and Bowhill 

(1965)] e 

3T 9 . ST 

n e kC v TT - Sin 1 A tK e - V L e (2) 


where K = the electron thermal conductivity 

e 

~3 

n^ ~ the electron density (cm ) 

G 

-i 

k = Boltzmann’s constant = 1.38x10 ergs/ K 

C v = the specific heat of the electron gas at a constant 

volume 

Q e = the electron heat production rate 

L " the electron heat loss rate 
e 

The electron heat conductivity is based on Spitzer’s formulation 
for a fully ionized Lorentz gas modified to account for the 
reduced electron mean free path due to the presence of the 

neutral gases. The main result of the modification is to reduce 

the effective thermal conductivity and increase the coupling 
to the neutral gas at lower altitudes. From Spitzer and H&°rm (1953) 
the thermal conductivity is given by 


K 


(s) 


= 1,23x10 


5/2 

x e 


ergs cm sec 


1 o 


K 


(3) 
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Introducing the modifications to account for the reduced mean 
free path 


K 


e 


K. 


TsJ 


+ 2 


1 .]-! 


n K 


n 


( 4 ) 


where K n is the thermal conductivity for electrons embedded 
in a proportionally dense neutral gas. Following Banks (1966) 
and Chapman and Cowling (1952) 

n_ 8 k T 

■>k t 
’n 


K e n " l%- )k [ 


n m. 


e 2 h 


(5) 


Q 


dn 


where Q dn is the velocity averaged momentum transfer cross section 

for electrons passing through a neutral gas. That is, 

™ 2 
— O r» 00 C HI V 

Q dn = [ '2 ' "'k T ] J v e q D (v e )exp(- ^®-)dv e (6) 


e 




N n = 


the velocity dependent momentum transfer cross 

section. This quantity is readily found from 

the published results of many experimenters 

(see McDaniel (1964) for a general reference) . 

t h 

the density of the n neutral species, 
the electron velocity. 


In the temperature range of interest, the thermal conductivity 
is given by 


A T 


5/2 


K e “ 


T 2T 
x e 

n 


5 


(ergs cm~ 1 sec~' 1 ' °K~"' 1 ') (7) 


S.. N Q -, 
m=l m dm 


1 + B 
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where 


A 

= 1 , 23x10”® 



B 

= 3.22xX0 4 



Q dl 

„1 <7 X 

= 2.82x10 T z -3. 41x10 

e 

-21 t 3/2 
e 

2 

cm 

Q d2 

= 2 .24x10 _16 +7.92x10" 18 ' 

F ^ 

"e 

„ 2 
cm 

Q d3 

■= 3.4xl0~ 16 


2 

cm 

Q d4 

= 5,6xl0 -16 


2 

cm 

Q d5 

= 5.47x10“ 15 -7.45x10" 19 

T 

e 

2 

cm 


The subscripts 1, 2 , 3, 4, 5 stand for the species Ng? 0^, 0, 

He, and H respectively. This ordering of the species is adopted 

throughout the remainder of the discussion. 

The rate of conversion of solar flux to thermal energy in 

the ambient electron gas, Q , is given by 

© 


Q = S N ± S s ± <4. a ^ 

e i=l 1 x lf X \ 


(I) 


exp(-T ) 


5 a 


(A) 


= S 


!,X 


X i=i cos x 


N. k T 

[ — 2. - C] 

m i B 


( 8 ) 

(9) 


where ^ = 

Cl) _ 


a 




the solar zenith angle 

the ionization cross section for species i and wavelength \ 


(A) 

a. J = the corresponding absorption cross section 
•*- ? A. 

0 = the solar flux in the range d\ about \ 

A 

s- = the energy conversion efficiency per electron-ion pair 

1 » A 

produced 
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g «=> the acceleration of gravity 
C a a small constant 

The values of the solar flux and cross sections are taken from 

Hinterreger et , al. (1965) in the wavelength range of 10-1027A. 

The energy conversion efficiency can be estimated from the 

available photon energy at each wavelength and the ionization 

potentials associated with each species. It is also necessary 

to consider the fraction of energy produced that is deposited 

locally in the electron gas by the photoelectrons. For the 

purposes of this study, an average value of 2 ev per electron- 

ion pair produced by photoionization has been found to be a 

suitable value. The effects of the detailed calculation to 

estimate e. will be presented in a future paper, Hanson and Cohen 
) A 

(1968) have attempted to experimentally determine such an efficiency 

factor from radar backscatter data for T and N . The result 

6 © 

of using such an efficiency is to shift the height of the maximum 
of both electron heat production and of the electron temperature 
upward in conformity with some of the observed data. 

The rate at which electrons lose their thermal energy to 
the ion and neutral gases, L Q , is determined from the velocity 
dependent momentum transfer cross section [Desloge (1962) ]. 

These cross sections, aside from the electron-ion coulomb cross 
section, are the result of laboratory measurements on nearly mono- 
energetic electrons passing through a background of the specified 
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neutral gas, The resulting velocity dependent cross sections 
are then averaged over a Maxwellian velocity distribution to 
obtain the energy loss from the passage of a hotter Maxwellian 
gas through a cooler one. Following Desloge (1962), the rate 
of energy loss from electrons to another species, s, is 


=4rm,_n 
es e s 


nhm 
e s 


-Qv. 


< ’ m e m s ) 5/,S ) (•» 5 ”“ v r 

2 ~7ZZrS727~ Z ~ Z \5 7a J v r <*o (v r )e 

' o 


dv. 


(m e +m s )<a (2nk) ^ m s^e +m e^s 


where n, m, and T refer to the appropriate densities, masses, 
and temperatures, and where 

v^ = the relative velocity between electrons and 

species, s 
T e T s 

n - [akC^ + j~) ] 

e s 

The experimental values of (v r ) can be approximated by 
a power series of the form 


n 


<l D (v r ) = “ A - v - + 2 


n n v r 


n 


n 


r Y on' 


( 11 ) 


The first term is used to fit the general trend of the data 
while the second can be used to account for important fluctuations 
a n is the strength of such a fluctuation at v p = v on > and is 
determined by the requirement that the total cross section, 
q^ , in a velocity range (A v) be given by 

f,v+Av 


v 


S A v n + £ a 8 (v -v ) ]dv 
£ n r £ n uv r on /J r 


( 10 ) 


The resulting form of the loss function is 


J es ~JtT 


T T i 
e . s - 2 


n 

IS r* i ,n+6- 


4kn n mm 

L -* — ^ 7rtf? [2k< =r * i ,r<T .- T s >t5 v " r<^> 

(m e +m s ) e s n 


+2 S a„v ® /2 exp(-n v„„) ] 


n 


n on 


on' 


( 12 ) 


For a multi constituent gas the electron heat loss function is 
assumed to be given by 


L = E + L . 
e es ei 


(13) 


Losses from electrons to ions via coulomb collisions are 

4 

In a 


n n e 

L . - 4a/2tf — k 
e i m m 4 


m m 

[k (- f + rp)] 3/2 
m A m . 

e i 


CVV 


(14) 


In A" 15.2 

From equations (12, (13) and (14) the loss terms may be expressed 
in the following form 


| L eg = 1.60xlO" 12 n e -f2.47xlO“ ;L8 [O](T e -T n )T e i 

+1 . 71xl0~ 19 [No ] (1-2 . 1 x10~ 4 T ) (T — T ) T 

4-1 c 6 XI © 

+1 . 21xl0~ 18 [0„ ] (1+3 . 6xlO -2 T (T -T )T ^ 

+3 . 1x 10“ 14 [N 2 ] (T e -T n )T“i+10 -13 [0 2 ] (T e -T n )T e - ^ 
+2.46xlO~ 17 [He]T e 2(T e ~T n ) 

—16 _ A i 

+9.63x10 [H](l~l. 35x10 T g ) (T e -T n )T 0 5 } (ergs cm' 

L e ^= 7.68x10 '*'^n e 2 (Tg-T^^) T^ 3 ^ 2 (ergs cro -8 see ” 4 ) 


ION TEMPERATURES 


The ion temperatures are obtained for each ionic species 
by solving the appropriate heat conduction equation in a manner 
exactly analogous to that for T . The source of thermal energy 
is provided mainly by the hotter electrons via coulomb collisions, 
while the losses arise from collisions with the neutral gas. 

This means that at low altitudes the temperatures for each 
ionic specie are strongly coupled to the neutral gas temperature, 
and at high altitudes to the electron temperature. Equality 
of the ion temperature, T i , at high altitudes with T e is prevented 
by the relatively large ion thermal conductivity which enables 
heat energy to be conducted downwards and ultimately dissipated 
by the neutral gas. The ion thermal conductivity, based 

on the ion-ion interaction, has been given by Chapman (1952) 
and has the same temperature dependence as K g . 

The ion temperature equation has the form 



St 


• 2 t . 
sxn I p 

n.kC *§0 L1 x 
x v 


ST 


Sz 


~] 


P.-L. 
x x 

n . kC 
x v 


(15) 


when the effects of motion are neglected. Since ionic composition 

. + 

is not being discussed in this work, only the major xon, 0 , 

is considered. However, it should be recognized that the usually 

minor ions, H + and He + , can be heated more efficiently than 0 + , 

and so can play an important role whenever their densities are 

+ + 

sufficiently close to the 0 density. Eor 0 , the following 
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expressions give K. , P. and L. . 

X X X 

K ± = 9.51xlO -20 T i 5/2 /ra i ^ (ergs cra“ 1 sec“ :L °K _1 ) 

-3 -1 

P^ = (ergs cm sec ) 

L, = n, (1.06 x10 -25 [N„]+9.3x10~ 26 [0„]+3.36x10 -27 [0](T,+T )® 

— —1 
+4.49x10 [He]+5.3xl0 [HlJOq-T ) (ergs cm sec ) 

NEUTRAL TEMPERATURE 

The heat conduction equation for the neutral temperature, 

T , is "based directly on that given by Harris and Priester (1962), 
but with three modifications. Since only the steady state results 
are of interest, the expansion velocity term is not included 
as it depends directly on the time rate of change of the neutral 
temperature. The heat input from the electron and ion gases 
are added to the heat input from the direct absorption of solar 
energy. This additional energy amounts to about 15-20% of the 
total heat input. Finally, a representative solar spectrum and 
corresponding ionization and absorption cross sections, as given 
by Hinterreger et . al. (1965), are used instead of the average 
values. Use of average values leads to a substantial error 
and can have a considerable effect on the distribution of thermal 
energy in the E and lower F regions of the ionosphere. 

The neutral gas temperature equation is given by 


r 
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BT 


n 


bt 


kN, 


T 


JL 

dz 


ST 

[K n aiT ] 


Q n~ L n ^ L es +L i 

T* 


kN, 


T 


kN, 


T 


(16) 


where 


K 


n 


SAN 

n mm A "i *. 

T n ergs cm‘“ 1 sec~ °K~ 


S N. 
m 


m 


is the density weighted approximation to the thermal conductivity 
of a multicomponent neutral gas 


A l= 18 ° 

A 2 =180 

A =180 
3 


A 4 =360 

A =2100 
5 

N =£ N c 
T m m vm 


c is the specific heat at constant volume 
vm ^ 

h 

for the m neutral constituent (c =3/2 

vm 

for atomic species and 5/2 for diatomic 
species) 

The energy input to the neutral gas from super elastic collisional 
conversion of the absorbed solar radiation is given by 

Q n - S ^ I x % (A > e X p(- Tx ) 

where B = the fraction of the absorbed energy that is transformed 
m\ 

into thermal energy. In this paper a constant value of 0.2 
is used so that realistic values of temperature could be 
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obtained for the range of solar flux considered in this 

paper. I is the intensity of the solar radiation in the 

X 

wavelength range AX about x • Finally, the loss term, L n , 
arising from the infrared radiation emitted by thermally excited 
atomic oxygen, is given by [Bates (1951), Lazerev (1963)]. 



[ 0 ] 


1.6x10’ 


18 


exp(-228/T n ) 


1+0 . 6exp (-228/T n ) +0 . 2exp (-325 . 3/T n ) 


/ “ 3 -lx 

(ergs cm sec ) 

(17) 


The distribution of the neutral constituents with altitude, z, 
is assumed to follow the hydrostatic distribution for z > z^. 


N(a) 

m 




T (a) 
n v ' 


exp {- 


a iTi g(z f ) 


m 


* kT(a') 
‘"'L n 


dz ’ ] 


(18) 


where the boundary values N (z T ) for the species [N 0 ], [0 o ], and 

m ju z z 

[0] are important input parameters whose variations are discussed in 
a later section. z T and z, are the lower and upper boundary 
altitudes, respectively. The remaining neutral species, He and 
H, are also described by the above equation. However, no 
significance can be attached to the values computed at low altitudes 
since it is only above about 400 km to 500 km that hydrostatic 
equilibrium is an appropriate approximation for He and H. Be- 
low this altitude they have little effect on the model because 
of their relatively small density. 
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ELECTRON DENSITY 


Once the various temperatures are specified by the preceding 
equations, the electron density can be obtained by solving 
the continuity equation together with the momentum equation. 

As with the ion temperatures , if the ionospheric composition 
is of interest, a separate continuity equation must be solved 
for each species, This involves specifying the expressions for 
the diffusion velocity and all important chemical and ionization 
processes. The solution of the resulting system of equations 
will be discussed in another paper. Since only the electron 
density is of interest here, a form of the continuity equation 
proposed by Shimizaki (1965) is used. This form accounts for 
the presence of the minor ions Ng*, 0 2 + ? and in E and 
lower F regions through modifications of the effective rate of 
photoionizat ion and chemical, recombination. The specific form 
used is 


a[o + ] 

~§t 


+ V.([0 + ]v 0+ ) 


s A. Q. - 

i i 

i=l 1 


L n 


(19) 


where 


Xf ^ 2 ]+X 2 ^2 ^ 


+ !i[0„] + A^[ No ] 


\2 


“I' 2 


a, 


A 4 *6. 


„e + (^ + ^ [0 ] + (^ + ^)[0 2 3 


Xo+Xc= X^(+Xo 

n e + (_L^)C0] + (-i-i)[fl 2 ] 
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and A 


j>Jl A a i— A - 1 

a 3 r 2 a x a 3 1 * 


The numerical values of the reaction rate coefficients are 
adopted from Whitten and Poppoff (1964)Ferguson (1968) , 
Warneck, (1966) , and are given in units of cm^sec ^ . 


a, “ 7xlO" 5 T “ 1 
1 n 


Xg = 2.5x10 


•10 


a„ = 1.5 xI 0~ 3 T ~ 3/2 
2 n 

x 4 

- 2xl0 -13 

0 

_ . ln -5 m “>3/4 
a Q " 5x10 T 
o n 

X 5 

= 2xl0~ 12 


X- « 2xl0 -11 

X 6 

= 2xl0~ 10 


X 2 = 2xl0~ 12 





These values were selected as representative from among 
the many choices in the literature. The variations between 
the various choices produce changes of little consequence for 
the F region electrons and 0 + ions, 

liquation (19) has been obtained from Shimizaki’s equation 
(13) (Shimizaki, 1965) by making several approximations. 

(From the numerical values of the rate coefficients it is 
possible to show that L/a^ , L/a 2 , and 6L/ a 3 become very small 
above about 200 km. Below 200 km the divergence terms, 

Vo (Co 2 + ]v 0 +) , Vo ([N 2 + ]v n +) , and V . ( [N0 + ]v N 0+ ) are almost 
entirely negligible. This means that the terms multiplied by 
L/cu , L/a 2 , and §L/a 3 need not be considered for steady state 
solutions . ) 
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Equation (19) does not include He* and H*. For the midday 
conditions discussed in this paper, the He* and H* densities 
are sufficiently less than the 0* density so that both the 
approximation [0*] « n and the neglect of He* and H* are 

V 

applicable throughout the height range of interest. Of 
course, this cannot be done when the diurnal variation is 
discussed since the influence of H* is usually dominant in 
the F region at night . 

The divergence term in equation (19) can be evaluated from 
the usual momentum equations [Burgers (1960) , Chandra (1964) ] . 


3v 

m e n e ( 3t- + *e>" 


“ ?P e +m e n ee- en e (E + c V B > 


- s K (v - V ) 
es e s 
s 


( 20 ) 


m 0+ [0 + ](-~ • + v Q+ -Vv 0+ ) = -7p 0+ +m 0+ [0" r ]g+e[0' r ](E + xB) 


t-i"? . _ r rt +' 


0 + 


s K 0 + s (v q+ - v s ) 


( 21 ) 


bv 


m nVaiT * v n* Vv n^ 


- yp +m N g - n (v - v ) 
*n n n & u ns v n s y 
s 


( 22 ) 


where 


P. = N.kT. 
3 3 3 


(23) 


g = the acceleration of gravity 


f 


e ** the electronic charge 
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B ” the earth's magnetic field in gauss 

K ij ~ K ji “ n ipij V ij 
Hid - ¥j /( “i +1 j) 

v. . = the momentum transfer collision 

frequency (see appendix for numerical 
evaluation) . 


It is now assumed that v.W can be neglected and that the 

neutral gas velocity v n can be treated as an unknown parameter 

of the problem combined with the hydrostatic approximation 

for the neutral gas distribution with altitude, Vp ^I^g. 

Finally, for the case where v = 0 and v . ./ (eB/m . ) «1 and 

^ ■* J J 

making use of the continuity equation V , (n e v e ~ [0 + ]v 0 -f )=0 , 
equation (22) for a horizontally stratified ionosphere may 


be approximated in the following form Chandra (1965) 


Bn 

e 

Bt~ 


■ 2 _ 
sm I 


A. 


■^[ n e k (T e +T . ) 3+ n e O e +m. ) g 

'-f J 0 + s v O + s + Hes v es ^ 


E aq 
n m nn 


L n 


(24) 


Equation (24) along with equations (2) , (15) , (16) , and 

(18) can be solved with the appropriate boundary conditions 

to give the profiles for N, n 0 , T 0 , T i , and T n in the height 

range between z T and z, . 

Xj u 


r 
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NUMERICAL METHOD 


The system of equations presented in the preceding sections 
is clearly not suited for finding an analytical solution under 
realistic conditions, Therefore, a particular numerical approach 
was adopted from among several possibilities that seemed feasible 
for this problem. Direct solution of the steady state equations 
was abandoned in favor of methods where the time derivative 
in each equation asymptotically approaches zero as time increases. 
The direct use of the steady state equations is very convenient 
when solving one of the equations of the system by itself. 

Rishbeth and Barron (1960) demonstrated in an ionospheric 
application how the difficulties in applying the boundary 
conditions can be overcome for the steady state equations by 
a trial and error procedure. However, when the method is 
applied to a system of coupled equations the result is far too 
cumbersome to be practical. The use of the time dependent 
equations to obtain an asymptotic approach to stationary 
conditions permits a straightforward numerical solution 
and leads naturally to the solution of the same equations for 
their diurnal behavior. 

All the equations are of the form 


a*i 


3 *i 


at?. it. 

+A i (a,t; + B jL (a Y ± ; t? +c(a,t; ; V „)-0 

i$0 


\ 
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where i and j are subscripts referring to the several dependent 

variables . 

* 

As far as the authors know, there is no theoretically 
described method for solving such a system of equations. 

i 

In the absence of such methods there are two basic practical 
approaches to be taken. First, one can linearize the equations 
so as to apply the well known theorems of linear and quasi- 
linear systems , [Ames (1965)2 . Second, one can adopt a partial 
linearization and iteration technique which is applied re- 
peatedly at each time step until the nonlinear equations are 
solved to some desired degree of accuracy. For example, the 
coefficients A, B, and C are calculated from an arbitrary set 
of initial values for they^. The solution is then obtained 
for the next time, 0 + At , by the implicit difference technique 
and iteration to obtain ^(O+At). ft ext , A, B, and C are re- 
calculated using the new values, and solutions are obtained 
for Y^(0+2At) . If the repeated process is non-divergent , stable, 
and yields a set of values, Y^(t) , that simultaneously satisfy 
the original system of equations and their boundary conditions 
to an accuracy within the truncation error caused by numerical 
differentiation, then a solution has been obtained. For 
the equations applied to the ionosphere, a great deal of care 
must be taken so that numerical significance is not lost when 
evaluating the space derivatives. This is especially true of 
the second derivatives in the height region from 100 to 200 km. 


i 
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The form of the boundary conditions applied are 


■= fj(t) 


sy. 

dsT 


= ^(t) 

u,t 


where 




e’ T „> T i’ n e> 


a L = 100 km 


« 1000 km 


f.(t) and 0 . (t ) are known specified functions 
1 1 

of time. 


Representative boundary conditions for the individual equations 
are shown in Table 1. 


Dependent Variable 


n 


TABLE 1 
Upper BC 


Lower BC 


velocity specified production=loss 


T e 

Sz - 0.5 /km 

T e 

= 226°K 

T i 

—• = 0 . 5°/km 

T i 

« 225. 5°K 

T n 

ST 

= 0 .0°/km 

T n 

= 225°K 


The choice of these values is discussed in detail in a later 


section . 
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DISCUSSION OF THE SOLUTIONS 

In this section a few selected cases are presented to show 

the general characteristics of the simultaneous solutions of 

the equations for T e , T i , T , and n e . All the cases presented 

are for conditions of latitude and dip prevailing at Wallops 

Is l and , Virginia (latitude of 37° and dip of 70°) or at Arecibo, 

Puerto Rico (latitude 18° and dip 50°) , and for a range of 

solar flux corresponding to high and low activity (noontime 

conditions with a declination competed for March 2) . 

To facilitate comparison between the various results, only 

the solar flux entering the ionosphere and the density of the 

neutral constituents at 100 km are varied as parameters of the 

problem. This means that the neutral atmospheric composition 

and corresponding heat and ionization input functions (Q . Q r , 

e ij 

Q , and P ) are determined in a self- cons isi ent manner with T , 
no © 

T • . T , and n_ and their boundary conditions . There are other 
free parameters that can be varied (reaction rates, collision 
cross sections, ionization cross sections, etc.), but to avoid 
an unnecessarily complex discussion they were held constant 
at values given by recent experimental results. Reference 
values for the N 9 , 0^, and 0 densities at 100 km were 
extrapolated from the measured values obtained 

during the Geoprobe# rocket flight of March 2, 1966 s [Hall et al 
1967] This flight made simultaneous measurements of the 

^Preliminary d~ata presented at the A.G.U. meeting in Washington, 
D„ Co (April 1968) by Brace, L. , Findlay, J„ A., and Mayr, 

H. G. (paper Ga89) and by Pelz, Do T. and Newton, G. P. 

(paper GAB 6) . 
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neutral composition, neutral temperature, ion composition, 
and electron density. Tho boundary values for He and H at 
100 km are given in Table 2. The gradient of 0.5°/km at an 
altitude of 1000 km for T Q and T^ was selected as being 
representative of the available data which ranges from 0°/km 
to 2°/km„ The boundary condition for the electron density 
was obtained by specifying the electron velocity as some 
reasonable value. It was found that velocities less than 
100 meters per second had no significant effect on the solutions, 
while those in excess of this value caused considerable dis- 
tortions from measured profiles. In addition, velocities 
of about 100 meters per second are an upper limit since they 
correspond to a flux of about 10 cm" sec . A flux of this 
size cannot be maintained for an extended period of time either 
into or out from the upper ionosphere[Geisler (1967)3 . The 
lower boundary at 100 km for the electron density was selected 
so as to satisfy the conditions of chemical equilibrium under 
all circumstances. If a different value was tried at 100 km, 
the solution at 110 km automatically assumed a value such that 
chemical equilibrium was reestablished. Significant deviations 
from chemical equilibrium began to appear at 120 km and higher. 
.Finally, the values of T , T., and T were selected as 226°, 
225.5°, and 225°K at 100 km, respectively. It was found that 
the solutions are relatively insensitive to this choice. There 
are three final parameters which help determine the character 


r 
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of the solutions. These are the neutral and charged particle 
heating efficiencies and the ionization efficiency. The first 
is an estimate of the fraction of solar energy absorbed by 
the neutral gas that becomes transformed into thermal (kinetic) 
energy by superelastic collisions involving neutral particles 
in an excited state. The second is an estimate of the number 
of electron-volts per ion-e rectron pair created that is 
ultimately deposited in the ambient electron gas as kinetic 
energy. This process involves the removal of energy from 
high velocity photoelectrons and their associated secondary 
impact ionization electrons. Several authors have made de- 
tailed calculations for this process [Dalgarno et „ al. (1963), 
Mariani (1964)]. Finally, the ionization efficiency accounts 
for the number of electrons produced per incident photon both 
by primary and secondary ionization process. The values adopted 
here are 20%, 2 ev/electron-ion pair, and 1.5 electron-ion 
pairs per primary electron ion pair produced. Actually, 
these quantities should be determined as self-consistent height 
dependent quantities along with the other unknowns upon which 
they depend. The uncertainty in the theoretical expressions 
for the efficiencies is too great to warrant the added complexity. 

Since the values measured by the Geoprobe rocket are used 
as a reference in obtaining the lower neutral atmospheric 
boundary conditions, a comparison is presented between the 
theoretical and experimental profiles. In Figure 2 the upper 
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half corresponds to smoothed curves drawn through the data 
points obtained from the Geoprobe rocket. The T n values are 
determined from the measured molecular nitrogen densities 
(the dashed portion of the curve is an extrapolation based 
on the assumed near isothermality at high altitudes) , the T 

© 

profile was taken from the cylindrical probe data, and the 
electron density from the propagation experiment carried on 
board the rocket. The lower half of the figure represents 
the corresponding theoretical model. As mentioned before, a 
fixed set of upper boundary conditions was used with no attempt 
tc match the measured gradient from the Geoprobe experiment 
of 2°/km 0 In spite of this, the agreement is quite reasonable 
over the entire altitude range for both T and n . Much closer 

II w 

agreement can be obtained when the fully time dependent 
equations are solved. The time dependent solutions admit the 
presence of large thermal gradients at high altitudes for 
limited portions of the day as long as the downward heat flow 
during the day does not exceed the available energy from the 
protonosphere. In the steady state models presented here the 
use of 2°/km as a boundary condition on at 1000 km leads 
to unrealistic results. The free parameters used to obtain 
the theoretical profiles are given in Table 2. The solar 
flux listed as 0.95*HF is scaled to be 95% of the flux given 
by Hinterreger, et . al. (1965). The adjustment of the solar 
flux is done so that the neutral exospheric temperature 
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TABLE 2 


Latitude = 37° 


Dip = 70 


o 


X = 39° (solar zenith 
angle) 

Day = March 2 


Solar Flux = 0.95*KF 


Neutral Density at 100 
[S 2 ]=4.65x10" 12 [0 2 ]=1.25xX0 12 

[He] =2. 5 Ox 10 8 


-3 

km (cm ) 
[0]=6.75xl0 l:l 

[H]=l . OOxlO 8 


approximately corresponds to the measured values. However, 
the resulting solutions showed that the adjusted solar activity 
necessary to give agreement with experimental results is close 
to the observed solar activity for March 2, 1966, It was also 
required that the range of solar flux used be approximately 
a factor of 2 in order to change the exospheric neutral 
temperature from solar minimum conditons (about 700°) to solar 
maximum conditions (about 1800°) . The necessary flux was 
obtained by multiplying the Hinterreger flux (HF) (corresponding 
to a 10.7 cm flux of about S=85) by a series of factors ranging 
from 0,85 to 1.75. 

Before discussing the effects of the varying solar flux, 
it is of interest to see how the neutral composition effects 
the solutions. In Figure 3 the upper half shows the results 
obtained by varying the molecular nitrogen density at 100 km 
for two values of the solar flux, 1,15*HF and 1,75*HF. The 
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main effect of this variation is to increase the neutral 
temperature and thereby the absolute concentration of atomic 
oxygen. Since this also implies an increase in optical depth 
and the race of ionization (and thus heat input) , it is con- 
sistent that the height of the electron temperature maximum 
rises as does the height of the maximum of electron density . 

An important change shows up in the gradients of and T i 
(T^ not shown) . These gradients control the shape of the 
electron density profile to a large degree. Of course, the 

profiles obtained for n at high solar flux cannot be realistic 

e 

because of the high altitude of the F2 peak [Rishbet h (1968)]. 

It is shown later that varying the atomic oxygen density 

at the lower boundary in proportion to the solar flux is 

sufficient to remove the difficulty. The lower boundary value 

11 -3 

used for atomic oxygen in figure 3 is 6.75x10 cm . The results 
and conclusions for varying the molecular oxygen density, 
as shown in the lower half of Figure 3, are analogous to those 
for molecular nitrogen. 

It should be pointed out that the effects of varying atomic 
oxygen depend on the way the variation is introduced. In the 
above case, the variation was caused by an increase in T n 
brought on by another source. However, if only the 

value of atomic oxygen at 100 km is increased, then a different 
sequence of events takes place. Specifically, the neutral 
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temperature decreases since the neutral heat loss function is 
proportional to the atomic oxygen density and the increased 
heat input from the electrons and ions is insufficient to 
make up the difference unless the increase in atomic oxygen 
density is small. In turn, the decrease in T n causes a 
depression of the peak heights in T 0 and n e , a decrease in 
n above the F2 maximum, and an increase below the maximum. 

Thus, more than any other neutral constituent, atomic oxygen 
plays the basic role in determining the characteristics of 
the E and F regions. This is particularly so, since it is also 
the neutral constituent that is most responsive to changes 
in the solar flux. Thus, variation of atomic oxygen with 
solar flux must play an important role in the F region behavior. 

Such effects should also be seen if there are other energy 

inputs into the ionosphere that significantly change the 

dissociation rates or raise the turbopause level (by expansion 

of the mesosphere due to heating). For example, the effects 

seen long after the onset of a magnetic storm may arise by the inverse of 

the atcve processes . This effect is discussed in the accompanying paper. 

The solar cycle variation of T , T . , T , and n at two 

€ jl n © 

latitudes was investigated for two sets of boundary conditions. 

The first set has the atomic oxygen density independent of the 
solar flux, while the second set assumes an empirical linear 
relationship between the two. Figures 4 to 9 show the effects 
of the solar flux variation (in steps of 0.85, 1,15, 1.45, and 


l 


r 
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1.75*HF) for three fixed values of the atomic oxygen density 
at 100 km, 5.39X10 11 , 6.75X10 11 , and 8 . llxlO l:L cm -3 . The first 
three of these figures are for a latitude and dip correspond- 
ing to Wallops Island, Virginia and the second three are for 
Arecibo, Puerto Rico. Some general features can be noted by 
direct comparison among the figures. First, in all cases the 
height of the F2 maximum increases with increasing solar flux, 
as does the neutral temperature, the total electron content, 
the total electron heat content (1.5 n kT dz) , the total 
ion and neutral heat contents, the height and magnitude of 
the maximum of T , and the magnitude of 'IV . These changes 
are accompanied by increases in the absolute concentrations of 
the neutral constituents at a fixed altitude, and by changes 
in the ratios of their densities. However, as can be seen 
from Figures 4 to 6 there are two distinct forms for the T 0 
profile. The first form shows a peak in the vicinity of 
200 km followed by a steep negative temperature gradient, while 
the second form is monotonically increasing except perhaps for 
a small region near 200 km. Of course, there exists a smooth 
transition between one form and the other, but over a moderately 
narrow range of the free parameters. As is shown, a change 
from 0.85 to 1.15*HF is sufficient to complete the transformation. 
The T profile labeled 1.05 in Figure 5 is included to show 
the smooth transition. The ion temperature associated with 
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the first form of T increases smoothly with solar flux in 

© 

the exosphere following the trends of T . At lower altitudes, 

below about 200 km, the ion temperature closely approaches T n 

in behavior and magnitude. However, as is shown in figures 

4 to 6, whenever the second form of T is encountered the 

© 

behavior of with solar flux is inverted until the solar 

flux 3b increased to where the first form of T appears. The 

effect of the second form of T 0 on the electron density profile 

is to strongly change the slope above the F2 maximum. The 

structure in the n Q profiles below the F2 maximum is caused 

by the gradients in the temperature profiles and by the various 

reaction rates contributing to the net electron production rate. 

The structure shown is suggestive of typically observed electron 

density profiles, but this comparison cannot be pursued until 

the detailed ion composition corresponding to this structure 

is also discussed. Finally, the effect of latitude shown in 

the two sets of figures is to cause the heights of the maximums 

of T and n_ to increase. That is, as the latitude decreases 
f i e 

the solar zenith angle decreases and the dip angle decreases, 
both of which have the effect of raising the heights of the 
temperature and density maximum. 

As mentioned before, the unrealistic-ally high altitudes 
c:f the F2 maximum for high solar flux values can be overcome 
by increasing the atomic oxygen content in proportion to the 
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solar flux. Such increases are consistent with those proposed 
by Shimizaki (1967) in a recent paper discussing this aspect 
of the neutral atmosphere. That is, the photodissociation rate 
of C>2 , eddy mixing, and the reactions into which atomic oxygen 
enters are such that the variation with solar flux is reasonable. 

The results of this procedure are shown in Figures 10 and 11. 

Thus, it is possible to obtain reasonable F2 maximum heights 
while still enabling the other parameters, T 0 , T^ , T , and 
n e , to change in a way appropriate to the solar flux. Because 
of the neutral temperature’s dependence on the atomic oxygen 
density for heat loss, the resulting T n values are lower in 
general for increased atomic oxygen density . The trend is 
not absolutely clear cut because of the increased ionization 
and thus the increased heat transfer from the electron and ion 
gases to the neutral gas. Even though the trends shown are 
consistent with observations, a detailed comparison is not 
made here since the atomic oxygen variation was carried out 
on a strictly empirical basis. The inclusion of the neutral 
chemistry and the associated reaction rates should make quantitative 
comparison possible. 

< 

CONCLUSION 

Theoretical profiles for T , T. , T , and n^ have been obtained 

e • i n e 

by numerical solutions of the coupled partial differential 
equations describing the E and F regions of the ionospheric 


l 
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plasma. The results obtained show the complex interrelationship 
between the several dependent variables and thus the necessity 
of obtaining simultaneous self-consistent solutions. In 
particular, the following conditions can be enumerated as the 
result of the present investigation, 

1. The neutral temperature of the upper atmosphere is 
strongly dependent on the solar flux and on the neutral 
composition at the turbopause level. The temperature increases 
with increases in the Og density or Ng density but decreases 
when the density of atomic oxygen is increased. 

2. The electron density above the F2 peak increases with solar 

flux as does h max , n max , and the total electron content. In addition, 
the electron density in the region of the F2 maximum increases with T 
unless the increase in T n is the result of a decrease in the 
atomic oxygen concentration at low altitudes, 

3. The composition at the lower boundary cannot be assumed 
constant throughout the solar cycle since it leads to unrealistic 
values of the F region parameters for sunspot maximum conditions. 

This difficulty was overcome by assuming an increase in the 
atomic oxygen density at the reference altitude proportional 

to the solar flux. The need for this variation shows the necessity 
for constructing a more complete model incorporating the neutral 
chemistry and mixing effect that are so important near the 
turbopause level. 


? 



33 


4, The electron temperature does not show a definite 

pattern of variation of the same regularity as the neutral 

temperature or the electron density. However, the integrated 

* 

heat content, c v n Q kT e dz, increases with solar flux. The 
electron temperature is characterized by two types of profiles. 
The first shows a peak in the vicinity of 200 km while the 
second type is monotonically increasing almost everywhere. 

The smooth but rapid transition between the two types of 
profiles when conditions are changed occurs during low solar 
activity. The transition is governed mainly by the value of 
T , going from a monotonically increasing profile to a peaked 
profile with a considerably reduced exospheric electron 
temperature. For high solar flux, T varies in proportion 
to the solar flux at all altitudes above the peak. 

5. The ion temperature behaves quite regularly unless 
the electron temperature undergoes a transition from one form 
to another as conditions are changed. At low altitudes the 
ion temperature closely follows the neutral temperature. 
However, above approximately 200 km, T^ starts to follow the 
trends in the electron temperature and frequently approaches 
T e at high altitudes . 
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A DDENDUM 

At the time of the completion of this work, Dalgarno 
and Degges (Planetary and Space Science , JU5, 125, 1968) 
published a paper showing that the dominant cooling mechanism 
fur the electron gas results from the fine structure excitation 
of atomic oxygen. We have since included this heat loss 
in our calculations. A typical example of the changes re- 
sulting from the additional electron cooling is shown in the 
following table. The conditions used correspond to the case 
HF-1.15 in Figure 5. 

Ionospheric Parameters Including Atomic Without Atomic Oxyg&r 

Oxygen Fine Fine Structure Cooling 

Structure Cooling 


T 

e 

(3000 km) 

°K 

1750 

1804 

T 0 

(500 km) 

°K 

1349 

1416 

T e 

(max) °K 


1195 

2274 

ll 

(T max) kra 


190 

180 

1 . 
'.V 

(1.000 km) 

°K 

859 

859 

T. 

(1000 km) 

°K 

1413 

1530 

n e 

(1000 km) 

-3 

cm 

3 . 2 6x10 3 

3 .98x10 

n e 

(500 km) 

-3 

cm 

6 . 47x10 4 

7.06x10 

n e 

(max) cm” 

3 

7.44xl0 5 

6.48x10 

h 

(n max) km 


260 

260 


As can be seen from the table, the main effect of this 
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cooling mechanism is to decrease T in the lower altitudes 

e 

while leaving all the other parameters substantially unchanged. 
The main conslusions of this paper and the accompanying paper 
("F-Region Ionization and Heating During Magnetic Storms") are 
not altered even though the electron temperature changes are 
considerable. The detailed implications of the atomic oxygen 
fine structure cooling of the electron gas will be discussed 
in a future paper. 
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APPENDIX 


The following are a list of effective collision frequencies 
for momentum transfer, v . . , taken from the equation for momentum 
transfer between two interpenetrating Maxwellian gases of 
possibly different temperatures and masses. The defining equation 
is 


v 




8 

37 fr 


„ 

/ij 




v q . . (v )exp(-Q. .v )dv 
r h gij r r ij r 7 


From this equation and the experimentally obtained velocity 
dependent momentum transfer cross sections, q di ^ , numerical 
values of v.. can be obtained (Banks, 1966) . 

A J 

For electron-neutral collisions, 

v el (e-N 2 ) =2 . 33xl0 -11 [Ng] f 1-1 . 21xlO -4 T e 3T e 
v e2 ( e -*° 2 )=l. 82 x 10 '' 10 [ 0 2][ 1 + 3 .6xlO -2 Te^]T e ^ 

v e3 (e->0)=2 .8xl0 -10 [0] {T 0 ^ } 

in A 

v e4 ( e ^ He )= 4 .6xl° [He]T e 2 

v e5 ( e -*H)= 4 . 5 xl 0_9 [ H ]{l-l. 35 x 10 “ 4T e 3T e 2 

For ion-neutral collisions, 

v~ (0 + -neutrals)=5.181xl0 -11 [0](T. +T ) ^ (1-0 .064 In (T. -IT )) 2 
on in in 

+1.798xl0~ 9 J N g (a g / (U 3 S >^ 

s 

s^3 

• ” i 
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o 

where the polarizabilities, a (in units of cm°) are given by 


-24 

0^=1. 76x10 ** 

a^=2 . 10x10*" 

■25 

-94 

a 2 =l. 60x10 

a c =6 . 70x10 
5 

•25 

a 3 =8.90xl0~ 25 
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